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1 Introduction 

In this chapter we review the precursors of catastrophic avalanches (global 
failures) in several failure models, namely (a) Fiber Bundle Model (FBM), 
(b) Random Fuse Model (RFM), (c) Sandpile Models and (d) Fractal Over- 
lap Model. The precursor parameters identified here essentially reflect the 
growing correlations within such systems as they approach their respective 
failure points. As we show, often they help us to predict the global failure 
points in advance. 

Needless to mention that the existence of any such precursors and detailed 
knowledge about their behavior for major catastrophic failures [1] like earth- 
quakes, landslides, mine/bridge collapses, would be of supreme value for our 
civilization. So far, we do not have any established set of models for these 
major failure phenomena. However, several reasonable models of failure have 
already been developed in various contexts. We review here some of the pre- 
cursors of global failures in these models. 
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2 Precursors in Failure models 

2.1 Composite material under stress: Fiber bundle model 



T 




Force/elongation 

Fig. 1: Fiber bundle model (left) and the force-response of a single fiber (right). 

Fiber bundle model represents various aspects of fracture-failure process 
in composite materials. A bundle of parallel fibers, clamped at both ends 
(Fig. 1 (left)) represents the model where fibers arc assumed to obey Hookean 
elasticity up to the breaking point (Fig. 1 (right)). The model needs three basic 
ingredients: (a) a discrete set of N elements (b) a probability distribution of 
the strength of fibers (c) a load-transfer rule. Peirce [2] initiated the model 
study in the context of testing the strength of cotton yarns. Since then this 
model has been studied from different views. Fiber bundles arc of two classes 
with respect to the time dependence of fiber strength: The 'static' bimdlcs 
[2-10] contain fibers whose strengths arc independent of time, where as the 
'dynamic' bundles [11-13] are assumed to have time dependent elements to 
capture the creep rupture and fatigue behaviors. According to the load sharing 
rule, fiber bundles are being classified into two groups: equal load-sharing 
(ELS) bundles and local load-sharing (LLS). In ELS bundles, intact fibers bear 
the applied load equally and in LLS bundles, the terminal load of the failed 
fiber is given equally to all the intact neighbors. With steadily increasing load, 
a fiber bundle approaches the failure point obeying a dynamics determined by 
the load sharing rule. The phase transition [7] and dynamic critical behavior of 
the fracture process in such democratic bundles has been established through 
recursive formulation [8, 9] of the failure dynamics. The exact solutions [9, 10] 
of the recursion relations suggest universal values of the exponents involved. 
Attempts have also been made [11] to study the ELS ans LLS bundles from 
a common framework introducing a single parameter which determines the 
load transfer rule. 

We discuss here, failure of static fiber bundles under global load sharing 
(ELS) for two different loading conditions: (a) load increment by equal amount 
and (b) quasi-static loading -which follows weakest link failure at each step of 
loading. We show analytically the variation of precursor parameters with the 
applied stress, which help to estimate the failure point accurately. 



Precursors and prediction of catastrophic avalanches 



3 



Recursive dynamics in ELS bundle: Precursors 

ELS model assumes that the intact fibers share the apphed load equally. The 
strength threshold of a fiber is determined by the stress value it can bear, and 
beyond which it fails. Usually, thrcsolds arc taken from a randomly distributed 
normalised density p(x) within the interval and 1 such that 

/ p{x)dx = l. (1) 
Jo 

The breaking dynamics starts when an initial stress a (load per fiber) is ap- 
plied on the bundle. The fibers having strength less than a fail instantly. Due 
to this rupture, total number of intact fibers decreases and effective stress in- 
creases and this compels some more fibers to break. Such stress redistribution 
and further breaking of fibers continue till an equilibrium is reached, where 
either the surviving fibers are strong enough to bear the applied load on the 
bundle or all fibers fail. 

This self organised breaking dynamics can be represented by recursion 
relations [8, 9] in discrete time steps. Lot Ut be the fraction of fibers in the 
initial bundle that survive after time step t, where time step indicates the 
number of occurrence of stress redistribution. Then the redistributed stress 
after t time step becomes 

and after t+1 time steps the surviving fraction of fiber is 

Ut+i = 1 - P{at); (3) 

where P{(Tt) is the cumulative probability of corresponding density distribu- 
tion p{x): P{(7t) = Jq* p{x)dx. Now using Eq. (2) and Eq. (3) we can write 
the recursion relations which show how at and Ut evolve in discrete time: 

and 

Ut+i = l-P{a/Uty,Uo = l. (5) 

At the equilibrium or steady state Ut+i = Ut = U* and at+i = <Tt = a* ■ 
This is a fixed point of the recursive dynamics. Eq. (4) and Eq. (5) can be 
solved at the fixed point for some particular distribution. At the fixed-point, 
the solutions [8, 9, 10] assume universal form 

a*{a)=C-{ac-af'^; (6) 
U*{a)=C+{a,-of'^- (7) 
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where <Jc is the critical stress and C is a constant. Prom the recursions and 
their sohitions wc can derive the following response quantities: 

Susceptibility: Which is defined as the amount of change in U when the 
external stress changes by a infinitesimal amount and can be derived as 

= \{<Jc-cT)-^;fi=\. (8) 

Relaxation time: This is the number of step the bundle needs to come 
to a stable state after the application of an external stress. Prom the solution 
of recursive dynamics we get [8, 9] near CTc 

rcx(ae-a)-«;0=^. (9) 



dU*{a) 
da 



Inclusive avalanche: This is the amount of avalanche per step of stress 
redistribution 

f = ^-^- (10) 

At a = the dynamics becomes "critically slow" as 

C/t~i"^;7 = i; (11) 



which suggests that at CTc inclusive avalanches follow a power law (exponent 
—2) decay with time. 



Prediction of global failure point 
(A) Using X aiid r 
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Fig. 2: Variation of and with applied stress for a bundle having N = 
50000 fibers and averaging over 1000 samples. We consider uniform distribution of 
fiber threshold. 
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We found that susceptibility (x) and relaxation time (r) follow power laws 
with external stress and both of them diverge at the critical stress. Therefore 
if we plot and with external stress, we expect a linear fit near critical 
point and the straight lines should touch X axis at the critical stress. We 
indeed found similar behavior (Fig. 2) in simulation experiments. 

For application, it is always important that such prediction can be done 
in a single sample. We have performed the simulation taking a single bundle 
having very large number of fibers and we observe similar response of x and 
T. The prediction of failure point is also quite satisfactory (Fig. 3). 
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Fig. 3: Variation of x and t ^ with applied stress for a single bundle having 
N = 10000000 fibers with uniform distribution of fiber threshold. 

(B) Using inclusive avalanche 





Fig. 4: Log-log plot of inclusive avalanche with step of load redistribution for 
7 different stress values [left] and log-normal plot of the same for 4 difi^erent stress 
values below the critical stress [right] . The simulation has been performed for a single 
bundle with N = 10000000 fibers having uniform distribution of fiber threshold. 

When we put a big load on a material body, sometimes it becomes im- 
portant to know whether that body can support the load or not. The similar 
question can be asked in FBM. We found that if we record the inclusive 
avalanche, i.e, the amount of failure in each load redistribution -then the pat- 
tern of inclusive burst clearly shows whether the bundle is going to fail or not. 
For any stress below the critical state, inclusive avalanche follow exponential 
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decay (Fig. 4 (right)) with time step and for stress values above critical stress 
it is a power law followed by a gradual rise (Fig. 4 (left)). Clearly at critical 
stress it follows a robust power law with exponent —2 that we already get 
analytically. 

(C) Using avalanche distribution 
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Fig. 5: Magnitude of avalanches with step of load increment in FBM. 

An "Avalanche" is the amount of failure occurs as the system moves from 
one stable state to the next stable state when the load is increased quasi- 
statically. We can simply measure it by counting the failure elements between 
two consecutive load increment. Such a series of avalanches has been shown 
in Fig. 5 up to the final failure point. If we record all the avalanches till fi- 
nal failure, the avalanche distribution follows a universal power law [4]. This 
avalanches can be recorded experimentally measuring the acoustic emissions 
during fracture-failure of materials. We want to study whether the avalanche 
distribution changes if we start gathering the avalanches from some interme- 
diate states of the breaking process (Fig. 6 (left)). In simulations we see that 
the exponent of avalanche distribution shows a crossover between two values 
(—5/2 ans —3/2) and the crossover point (length) depends on the starting 
position (xq) of our measurement (Fig. 6 (right)). 
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Fig. 6: Average force versus elongation curve (left) and simulation results of 
avalanche distributions for different starting position (right) in FBM. 
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Now we are going to explain the above observation analytically for quasi- 
static loading. For a bundle of many fibers the average number of bursts of 
magnitude A is given by [4] 



D{A) A 



N 



A\ 



p{x) [1 - xp{x)/Q{x)] [xp{x)/Q{x)] 



A-l 



(12) 



X exp[-Axp{x)/Q{x)]dx, (13) 

where Q{x) = p{y) dy is the fraction of total fibers with strength exceeding 
X and Xc is the critical value beyond which the bundle fails instantly. For uni- 
form distribution (having upper boimd Xm), if we start recording avalanches 
from a intermediate point xo, the avalanche distribution becomes 



D{A) _ 
N ~ A\{xm-xo) 



-2x 



dx. (14) 



Introducing the parameter e = ^'L ^° and a new integration variable z 



D{A)_ 

N 



, we obtain 



A\{l + 2e) 



4/(l+2c 



(1 - ez)(2- ez)' 



For small e, i.e., close to the critical threshold distribution, we can expand 



ez + ln(l - ez) = -\e^z^ - \e^z^ + . 



(16) 



with the result 



D{A) _ A^-^e-^e^ „ _^,v/2 . Zi^-^-^ 



N 2A\ 



e-^^'^'l- zdz = (1 - e-«^^^) . (17) 



Using Stirling approximation A\ ~ A'^e "^V^ttA, this becomes 
D{A 



~ (87r)-V2z\-5/2 (1 - e-^/^=) , (18) 



N 

with 

1 _ x^ 
8e^ 8{xc — Xo 

Clearly, there is a crossover [5] at a burst length around Ac, so that 



^c = ^= .,„ (19) 



D{A) f (8/7r)V2e2 {ovA^A^ 



r (8/7r)V^ 
I (Stt)-^ 



N ~ I (87r)-VM-5/2 for Zi > Zi, 



(20) 



For Xo = OAOxrrn we have Ac — 12.5 uniform distribution. Wc found a clear 
crossover near A = Ac = 12.5 in the simulation experiment (Fig. 7). 
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Fig. 7: Crossover in avalanche distributions: bundle with 10'' fiber and Xq = 0.8xc 
(averaging over 20000 samples) [left] and a single bundle with 10^ fiber; xq = 
(squares) and xo = 0.9xc (circles) [right]. Dotted straight lines are the best fits to 
the power laws. 

The simulation results shown in the figures arc based on averaging over 
a large number of samples. For applications it is important that crossover 
signals are seen also in a single sample. Wc shovif (Fig. 7) that equally clear 
power laws arc seen in a single fiber bundle when is large. 

This crossover phenomenon is not limited to the uniform threshold distri- 
bution. The ^ = 3/2 power law [5] in the burst size distribution dominates over 
the ^ = 5/2 power law whenever a threshold distribution is non-critical, but 
close to criticality. Therefore, the magnitude of the crossover length correctly 
inform us how far the system is from the global failure point. 

2.2 Electrical networks within a voltage difference: Random fuse 
model 




Fig. 8: Random fuse model (left) and the current-response of a single Ohmic 
resistor (right). 
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Random fuse model [1] describes the breakdown phenomena in electrical 
networks. It consists of a lattice in which each bond is a fuse, i.e., an Ohmic 
resistor as long as the electric current it carries is below a threshold value. If 
the threshold is passed, the fuse burns out irreversibly (Fig. 8 (right)). The 
threshold (th) of each bond is drawn from an uncorrelated distribution p{th). 
All fuses have the same resistance. The lattice is a two-dimensional square one 
placed at 45° with regards to the bus bars (Fig. 8 (left)) and an increasing 
current is passed through it. Numerically, the Kirchhoff equations are solved 
at each point of the system. 
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Fig. 9: Avalanches in fuse model (left) and the ultimate fractured state (right). 

For the entire failure process avalanches of different sizes occur (Fig. 9 
(left)) and the system finally comes to a point where no current passes through 
the system (Fig. 9 (right))- we call it the global failure point. The average 
avalanche size has been measured [7] for such system and it follows a power 
law with increasing current or voltage: 

m ~ (/c - or TO - (V; - V)-^';-/ = 1/2 

Therefore if we plot with J or we can expect a linear fit which 
touches X axis at the critical current value Ic or critical voltage value Vc 
(Fig. 10). 
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Fig. 10: Plot of TO ^ with V (left) and the crossover in avalanche distributions 
in a fuse model (right). 
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If we record all the avalanches, the avalanche distribution follows a univer- 
sal power law with ^ ~ 3. But what will happen if we start recording the burst 
at some intermediate state of the failure process? We observe [5] that for a 
system of size 100 x 100, 2097 fuses blow on the average before catastrophic 
failure sets in. When measuring the burst distribution only after the first 2090 
fuses have blown, a different power law is found, this time with ^ = 2. After 
1000 blown fuses, on the other hand, ^ remains the same as for the histogram 
recording the entire failure process (Fig. 10). 

2.3 SOC models of sandpile 
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Fig. 11: A real sandpile and the sandpile model [16] on a square lattice. 

Growth of a natural sandpile is a nice example of self-organised critical- 
ity [15-19]. If sand grains are added continuously on a small pile, the system 
gradually approaches towards a state at the boundary between stable and un- 
stable states where the system shows long-range spatio-temporal fluctuations 
similar to those observed in equilibrium critical phenomena. This special state 
has been identified as the critical state of the pile where the response to ad- 
dition of sand grains becomes unpredictable: Avalanche of any size is equally 
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probable at this state. Therefore, we can expect system spanning avalanche 
(global failure) only at this critical state. 

BTW model and Manna model 

Bak et. al. [15] proposed a sandpilc model on square lattice which captures 
correctly the properties of a natural sandpile. At each lattice site there 
is an integer variable hij which represents the height of the sand column at 
that site. A unit of height (one sand grain) is added at a randomly chosen 
site at each time step and the system evolves in discrete time. The dynamics 
starts as soon as any site has got a height equal to the threshold value 
{hih= 4): that site topples, i.e., hij becomes zero there, and the heights of 
the four neighbouring sites increase by one unit 

hi,3 hij — 4, hi±ij hi±ij + 1, hij±i — » hij±i + 1. (21) 

If, due to this toppling at site any neighbouring site become unstable 

(its height reaches the threshold value), the same dynamics follows. Thus 
the process continues till all sites become stable {hij < hth for all 
When toppling occurs at the boundary of the lattice, extra heights get off the 
lattice and are removed from the system. With continuous addition of unit 
height (sand grain) at random sites of the lattice, the avalanches (toppling) 
get correlated over longer and longer ranges and the average height (hav) of 
the system grows with time. Gradually the correlation length (^) becomes of 
the order the system size L as the system attains the critical average height 
hc{L). On average, the additional height units start leaving the system and 
the average height remains stable there (see Fig. 12 (a)). The distributions of 
the avalanche sizes and the corresponding life times follow robust power laws 
[17], hence the system becomes critical here. 

We can perform a finite size scaling fit hc{L) = hc{oo) + GL~^/'^ (by 
setting ^ ~ I hc{L) — /ic(oo) \~"= L), where C is a constant, with v ~ 1.0 
gives he = hc{oo) ~ 2.124 (see inset of Fig. 12 (a)). Similar finite size scaling 
fit with I' = 1.0 gave hc{oo) ~ 2.124 in earlier large scale simulations [16]. 
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Fig. 12: Growth of BTW (a) and Manna (b) sandpiles. Inset shows the system 
size dependence of the critical height. 

Manna proposed the stochastic sand-pile model [16] by introducing ran- 
domness in the dynamics of sand-pile growth in two-dimensions. Here, the 
critical height is 2. Therefore at each toppling, two rejected grains choose their 
host among the four available neighbours randomly with equal probability. Af- 
ter constant adding of sand grains, the system ultimately settles at a critical 
state having height he- A similar finite size scaling fit hc{L) = hc{oo) + CL~^/'' 
gives ly ~ 1.0 and he = hc{oo) ~ 0.716 (see inset of Fig. 12 (b)). This is close 
to an earlier estimate he — 0.71695 [19], made in a somewhat different version 
of the model. The avalanche size distribution has got power laws similar to the 
BTW model, however the exponent seems to be different [16, 17], compared 
to that of BTW model. 

Sub-critical response: Precursors 

We arc going to investigate the behavior of BTW and Manna sandpiles when 
thay arc away from the critical state {h < /ic). At an average height hav, when 
all sites of the system have become stable (dynamics have stopped), a fixed 
number of height units hp (pulse of sand grains) is added at any central point 
of the system [8, 18]. Just after this addition, the local dynamics starts and it 
takes a finite time or iterations to return back to the stable state {hij < hth for 
all (i, j)) after several toppling events. We measure the response parameters: 
A — > number of toppling t number of iteration and ^ — > correlation length 
which is the distance of the furthest toppled site from the site where hp has 
been dropped. 
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(A) In BTW model 



We choose h„ 



4 for BTW model to ensure toppling at the target 



site. We found that all the response parameters follow power law as he is 
approached:Z\ oc {he — hav)^^, t oc {he — hav)^^, ^ gc {he — hav)~''', A = 2.0, 
II ^ 1.2 and ly ^ 1.0. Now if we plot A"'^^^, r"!/^ and ^''^/•^ with hav all 
the curve follow straight line and they should touch the x axis at hav = he- 
Therefore by a proper extrapolation we can estimate the value of he and we 
found he = 2.13 ± .01 (Fig. 13) which agree well with direct estimates of the 
same. 
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Fig. 13: Precursor parameters and prediction of critical point in BTW model. 
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(B) In Manna model 



Obviously we have to choose hp — 2 for Manna model to ensure toppling at 
the target site. Then wc measured all the response parameters and they seem 
to follow power law as he is approached:Z\ oc {h,. — hav)^^, t oc {he — hav)~'^, 
^ oc (/ic - havy-, A = 2.0, ^l ^ 1.2 and v ^ 1.0. As in BTW model, aU the 
curve follow straight line if we plot Z\~^/'^, r^^/^ and ^^^/"^ with hav and 
proper extrapolations estimate the value of he = 0.72 ± .01 which is again a 
good estimate (Fig. 14). 
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Fig. 14: Precursor parameters and prediction of critical point in Manna model. 



Our simulation results suggest that although BTW and Manna models 
belong to different universality class with respect to their properties at the 
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critical state, both the models show similar sub-critical response or precursors. 
A proper extrapolation method can estimate the respective critical heights of 
the models quite accurately. 

2.4 Fractal overlap model of earthquake 

It has been claimed recently that since the fractured surfaces have got 
well-characterized sclf-affine properties, the distribution of the clastic en- 
ergies released during the slips between two fractal surfaces (earthquake 
events) may follow the overlap distribution of two self-similar fractal sur- 
faces [20, 21, 22, 23]. To support this idea, Chakrabarti and Stinchcombe [24] 
have analytically shown using renormalization group technique that for reg- 
ular fractal overlap (Cantor sets and carpets) the contact area distribution 
follows power law. This claim has also been verified by extensive numerical 
simulations [25]. If one cantor set moves uniformly over other, the overlap 
between the two fractals (Fig. 15) change quasi-randomly with time. In this 
section we analyse the time series data of such overlaps to find the prediction 
possibility of a next large overlap. 

(a) (b) 



n=l 

rz IT 

n=3 n=3 

Fig. 15: (a) A regular Cantor set of dimension In 2/ In 3; only three finite gener- 
ations are shown, (b) The overlap of two identical (regular) Cantor sets, at n = 3, 
when one slips over other; the overlap sets are indicated within the vertical lines, 
where periodic boundary condition has been used. 

The time series data analysis 

We consider now the time series obtained by counting the overlaps m{t) as a 
function of time as one Cantor set moves over the other (periodic boundary 
condition is assumed) with uniform velocity. The time series are shown in Fig. 
16., for finite generations of Cantor sets of dimensions In2/ln3 and In4/ln5 
respectively. 

We calculate the cumulative overlap size Q{t) = mdt. We see that 'on 
average' Q{t) is seen to grow linearly with time t for regular Cantor sets. 
This gives a clue that instead of looking at the individual overlaps m{t) series 
one may look for the cumulative quantity. In fact, for the regular Cantor set 
of dimension In 2/ In 3, the overlap m is always 2*^, where k is an integer. 
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However the cumulative Q{t) = X]*^q 2^^* can not be easily expressed as any 
simple function of t. Still, we observe Q{t) ~ ct, where c is dependent on M. 
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Fig. 16: The time (t) series data of overlap size (m) for regular Cantor sets: (a) of 
dimension In 2/ In 3, at 8th generation: (b) of dimension In 4/ In 5, at 6th generation. 

We first identify the 'large events' occurring at time ti in the m{t) series, 
where m{ti) > M, a pre-assigned number. Then, we look for the cumulative 
overlap size Q{t) — /j*'^^ mdt, where the successive large events occur at times 
ti and ti+i. The behavior of Qi with time is shown in Fig. 17 for regular cantor 
sets with df =ln2/ In 3 at generation n — 8. Similar results are also given for 
Cantor sets with df =ln4/ In 5 at generation n = 6 in Fig. 18. ft appears that 
there are discrete (quantum) values up to which Qi grows with time t and 
then drops to zero value. 





Fig. 17: The cumulative overlap size variation with time (for regular Cantor sets 
of dimension In 2/ In 3, at 8th generation), where the cumulative overlap has been 
reset to value after every big event (of overlap size > M where M = 128 and 32 
respectively). 
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Fig. 18: The cumulative overlap size variation with time (for regular Cantor sets 
of dimension In4/ln5, at 6th generation), where the cumulative overlap has been 
reset to value after every big event (of overlap size > M where M = 2400 and 
2048 respectively). 

Finally we found that if one fixes a magnitude M of the overlap sizes m, 
so that overlaps with m > M arc called 'events' (or earthquake), then the 
cumulative overlap Qi grows linearly with time up to some discrete quanta 
Qi = lQQ, where Qo is the minimal overlap quantum, dependent on M and I 
is an integer. 



3 Conclusions 

Knowledge of precursors sometimes help to estimate precisely the location of 
the global failure or critical point through a proper extrapolation procedure. 
Therefore precursors which are available long before the global failure, can be 
used to resist an imminent global failure. 

In all the dynamical systems studied here, we find that long before the 
occurrence of global failure, the growing correlations in the dynamics of con- 
stituent elements manifest themselves as various precursors. In Fiber Bundle 
Model, the breakdown susceptibility x and the relaxation time r, both diverge 
as the external load or stress approaches the global failure point or critical 
point. The distribution of avalanches exhibits a crossover in power law ex- 
ponent values when the system comes closer to the failure point. Also the 
pattern of inclusive avalanches can tell us whether a bundle can support the 
applied load or not. In Fuse model, divergence of susceptibility and crossover 
in avalanche power law exponent are the signature of imminent breakdown of 
the system. Whereas in both BTW and Manna sandpile models the number 
of toppling A, relaxation time t and the correlation length ^ grow and diverge 
following power laws as the systems approach their respective critical points 
he from the sub-critical states. Therefore these parameters directly help to 
predict the critical point in advance. However in fractal overlap model the 
time series data analysis suggest that the cumulative overlap sizes can assume 
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some quantized values which indirectly helps to speculate whether a large 
overlap (event) is imminent or not. 
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